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Dielectric resonators are open systems particularly interesting due to their wide range of appli- 
cations in optics and photonics. In a recent paper [PRE 78, 056202 (2008)] the trace formula for 
both the smooth and the oscillating parts of the resonance density was proposed and checked for the 
circular cavity The present paper deals with numerous shapes which would be integrable (square, 
rectangle, and ellipse), pseudo-integrable (pentagon) and chaotic (stadium), if the cavities were 
closed (billiard case). A good agreement is found between the theoretical predictions, the numerical 
simulations, and experiments based on organic micro-lasers. 
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I. INTRODUCTION 

Open quantum (or wave) systems are rarely integrable 
and therefore difficult to deal with. Over recent years, 
this field of research has raised many crucial questions 
and various systems have been investigated. Here we 
consider open dielectric resonators for their wide range 
of applications in optics and photonics p], Q . In a first 
paper [3], the trace formula for these systems was de- 
rived in the semi-classical regime to infer their spectral 
features. More specifically in that paper both the expres- 
sions for the weighting coefficients of the periodic orbits 
and the counting function N(k) (mean number of reso- 
nances with a real part of the wave number less than k) 
were obtained and demonstrated analytically for two in- 
tegrable cases, the two-dimensional (2D) circular cavity 
and the ID Fabry-Perot resonator. In the present paper, 
we consider in detail 2D dielectric cavities with different 
shapes where no explicit exact solution is known. We 
compare the predictions of formulae obtained in [3| with 
numerical simulations and experiments based on organic 
micro-lasers. 

Resonance problems can be seen as counterparts of the 
scattering of an electromagnetic wave on a finite obstacle. 
This point of view turns out to be particularly interest- 
ing since such scattering problems have been extensively 
studied (see e.g. [4]). Rigorous results for the scattering 
of a wave on convex obstacles with Dirichlet boundary 
conditions were proved in [5] . Physical approach to these 
problems has been discussed in J6j . More recently some 
theorems were demonstrated in [JQ. The general struc- 
ture of the resonance spectrum on a transparent smooth 
obstacle was studied in Q. 

This paper is focused on careful investigations of spec- 



tral properties for 2D convex dielectric resonators, which 
are the open counterparts of the so-called 'quantum bil- 
liards'. The outline of the paper is the following. The 
formulas obtained in [3] are recalled and the numerical 
and experimental techniques are described in Sec. [TTJ 
Then different cavity shapes are explored and their prop- 
erties are compared with what is known for billiards. The 
square, rectangle, and ellipse cases are gathered in Sec. 
Hill We call such shapes 'regular shapes' since the cor- 
responding billiard problems are separable. In Sec. HVj 
the pentagonal dielectric cavity was chosen to illustrate 
a pseudo-integrable system. Eventually in Sec. [Vj the 
Bunimovich stadium is investigated as an archetype of 
a chaotic system. For completeness in Appendix \K\ the 
derivation of the Weyl's law is briefly presented. 



II. BACKGROUND: THEORY, NUMERICS, 
AND EXPERIMENTS 



Real dielectric resonators are three-dimensional (3D) 
cavities requiring that the 3D vectorial Maxwell equa- 
tions are used. When the cavity thickness is of the 
order of the wavelength, this problem can be approxi- 
mated to a 2D scalar equation following the effective in- 
dex model, which is widely used in photonics (see e.g. [l| 
and references therein). This approach has been proved 
to be quite efficient for our organic micro- lasers llOL 11 . 
Briefly, it assumes that the electromagnetic field can be 
separated in two independent polarizations, called TM 
(resp. TE) if the magnetic (resp. electric) field lies in 
the plane of the cavity (xy) [33[. In this 2D approxima- 
tion the Maxwell equations are reduced to the Helmholtz 
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equation 



(A xy + n 2 k 2 )4> = 



(1) 



where i/j stands for the z-component of the electric (resp. 
magnetic) field in TM (resp. TE) polarization. After res- 
olution, all the components of the electromagnetic fields 
can be inferred from ip. In Eq. (pQ) k is the wave number 
and n the effective refractive index. It is worth high- 
lighting that the error of this approximation is not well 
controlled [l2|. 

The boundary conditions in this 2D approximation are 
the following: 
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where v is a direction normal to the boundary and ipi 
(resp. ^2) corresponds to the field inside (resp. outside) 
the cavity. In the case of an open system such as a dielec- 
tric cavity, the resonances are defined as the solutions of 
(PQ) with the outgoing boundary condition at infinity: 
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Then the resonance eigenvalues, /c^, are complex with 
negative imaginary part 



(3) 



E is called the energy of the resonance whereas r is its 
lifetime. The wave numbers of the low- loss resonances 
(higher quality factors) are thus located closed to the 
real axis. 



A. Semiclassical trace formula 

Here for simplicity, we consider only TM polarization 
where the functions and their normal derivative are con- 
tinuous on the cavity boundary. In this case, it appears 
that the resonance spectrum splits into two subsets, de- 
pending on the imaginary part of the wave numbers. For 
one of the subsets, the wave numbers lie above a bound- 
ary 



7n 



< Im k n < 



(4) 



where ^/max is a certain constant which depends on the 
cavity, and the corresponding wave functions are mainly 
concentrated inside the cavity. These resonances are 
similar to the so-called Feschbach resonances. For the 
second class of resonances (called shape resonances) the 
wave functions are mainly supported outside the cavity 
and the corresponding eigenvalues have large imaginary 
parts. For smooth convex obstacles, it was shown (see 



e.g. |13| 
equality 



and references therein) that they obey the in- 



Im k n < — Const |Re k n 
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Hereafter we will focus only on Feschbach (inner) reso- 
nances, since they are the most relevant for lasers and 
photonics applications. They will simply be referred to 
as "resonances" from now on. 

The spectral density can formally be separated in two 
contributions 



d(k) =d(k)+d iosc \k) . 



(5) 



d(k) stands for the smooth part and is usually written 
through the counting function d(k) = dN(k)/dk which 
counts how many resonances in average have a real part 
less than k [34[. The oscillating part, d^ osc \ can be re- 
lated in the semiclassical regime kl ^> 1 (I is any charac- 
teristic length of the cavity) to a sum over the classical 
periodic orbits [14j . 

In [3], the semiclassical trace formula for open dielec- 
tric cavities was derived. It states that the counting func- 
tion of dielectric resonators can be written as follows: 



N(k) = n [ 
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where A is the area of the cavity, £ its perimeter, and 
f(n) a function of the refractive index involving elliptic 
integrals: 



f(n) = 1 + 



with 
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The derivation of (j6j) and details on r(n) and R(i) are 
given in Appendix [A] In the following, we will compare 
for various shapes the prediction of (j6]) to the function 
N(k) inferred from numerical simulations, and show a 
good agreement in all considered cases. In particular, we 
will stress the non trivial linear coefficient 



th 



r(n)£/47r 



(7) 



In this paper, in general n = 1.5, and so f (1.5) = 1.025. 

The oscillating part of the trace formula is written as 
a sum over classical periodic orbits (p.o.): 



d(° sc \k) = J2( c P einklp + c - c -) 



(8) 



p.o. 



where l p is the length of the orbit and c p its amplitude 
which depends only on classical quantities. We count a 
periodic orbit and the corresponding time-reversed or- 
bit as a single orbit. The expressions for the c p can be 
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derived in a standard way (see eg. |l4ll5j|). using the for- 
mula of the reflected Green function given in Appendix 
lAl As for billiard, it depends whether the orbit is isolated 
(i.e. unstable) or not. For an isolated periodic orbit 
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where M p , and R p are respectively the monodromy 
matrix, the Maslov index of the orbit, and the product of 
the Fresnel reflection coefficients at all reflection points. 
For a ray with an angle of incidence Xi the TM Fres- 
nel reflection coefficient at a planar dielectric interface 
between a medium with a refractive index n and air is: 



n cos x 



Rtm(x) = 

n cos x 

For a periodic orbit family 
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where A v is the area covered by the orbit family and (R p ) 
stands for the average of the Fresnel reflection coefficient 
over the family. 

Hereafter, to compare these theoretical predictions to 
numerical simulations, we will rather consider the Fourier 
transform of the spectral density d(l) in order to reveal 
the oscillating part: 



d(0 = E' 



-ik n l 



(12) 



where the k n are the complex eigenvalues calculated from 
numerical simulations. This function can be obtained 
from experiments as well. 



B. Numerical simulations 

The numerical simulations are based on the Bound- 
ary Element Methods which consists of writing the so- 
lution of (pQ) as integral equations on the inner and 
outer sides of the boundary and of matching them using 
the boundary conditions. The complex spectrum and 
the resonance wave functions (sometimes called quasi- 
stationnary states or quasi-bound states) are inferred 
from the obtained boundary integral equation. In ac- 
cordance with the experiments presented here based on 
polymer cavities, we used n = 1.50 inside the cavity and 
n = 1 outside (air). 



C. Experiments 

Dielectric resonators are widely used in photonics for 
fundamental research [16] and practical applications Q. 




FIG. 1: (Color online) Photographs through an optical micro- 
scope of some micro-lasers used in this paper. The in-plane 
scale is about 100 /im. 



Furthermore their wavelength range is not limited to op- 
tics and cover other electromagnetic domains like mi- 
crowaves [13] or Tera-Hertz waves [18j. Here we con- 
sider quasi-2D organic micro-lasers since they proved to 
be quite efficient to test trace formulae [l(| El • 

The cavities were etched by electron-beam lithogra- 
phy into a polymethylmethacrylate (PMMA) layer doped 
with a laser dye [35] spincasted on a silica on silicon 
wafer. This technology offers appreciable versatility in 
terms of shapes, while ensuring a small roughness and a 
good quality for the corners with a resolution better than 
a tenth of wavelength (see Fig. 6.5 in [11]). Some pho- 
tographs of the cavities studied in this paper are shown 
in Fig. [TJ The fabrication process is relatively fast and 
reproducible. At the end, the cavity thickness is about 
600 nm, while the in-plane scale is of the order of a few 
dozens of microns, which allows to apply the effective 
index approximation and therefore to consider these cav- 
ities as 2D resonators. 

The chosen cavity was uniformly pumped from above 
at room temperature and atmosphere with a frequency- 
double pulsed Nd:YAG laser (30 or 700 ps) and its emis- 
sion, integrated over 30 pump pulses, was collected side- 
ways in its plane with a spectrometer (Acton SpectraPro 
2500i) coupled to a cooled CCD camera (PixislOO Prince- 
ton Instruments). The spectral range of the emitted light 
depends on the dye laser. Here, for DCM, it is centered 
around 600 nm and so the kl parameter varies from 500 
to 1000. Consequently this experimental system is work- 
ing far away within the semi-classical regime while its 
coherence properties are ensured by lasing. 

The laser emission is mostly TE polarised [l9| , but for 
the features which are compared here with theory and nu- 
merics, there is not any predicted difference between the 
TE and TM cases. Among the resonator shapes studied 
in this paper, the pump polarisation plays a prominent 
role only for the square where it will be further developed. 
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For the other shapes, it will not be mentioned. 

As this paper is focusing on spectral features, we will 
consider only the emission spectra which, by default, 
were registered in the direction of maximal emission (i.e. 
parallel to the sides for square and pentagon [20| , parallel 
to the shortest axis for the rectangle [20|, and at an an- 
gle depending on the shape parameter for stadiums |2l|). 
Moreover in order to be closed enough to the theoretical 
case of a passive resonator, the cavities were pumped just 
above the laser threshold. Mode (and orbit) competition 
is then reduced. 

The typical laser spectrum is made of one or several 
combs of peaks connected (in a crude approximation) 
with certain periodic orbits. As shown in [10] the geo- 
metrical lengths of the underlying periodic orbits can be 
inferred from the Fourier transform of the experimental 
spectrum, which is an equivalent of the length density 
d(l). For instance, for a Fabry-Perot resonator of width 
a, the geometrical length of the single periodic orbit is 
L = 2a and the dephasing after a loop should be a mul- 
tiple of 2tt: knL = 2i\m with m G N*. Then the spacing 
between the comb peaks verifies Ak = 27r/nL, leading to 
a periodic comb pattern and a Fourier transform of the 
spectrum peaking at nL. With our experimental set-up, 
the precision on the geometrical length L reaches 3 % 
after taking duly account of the dispersion due to the ef- 
fective index and the absorption of the laser dye. So the 
refractive index which should be used to interpret the 
Fourier transform is 1.64 for these actual experiments (it 
is different from the bulk refractive index 1.54 and the 
effective refractive index 1.50) [lOj. 



III. REGULAR SHAPES 

This section deals with square, rectangle, and ellip- 
tic dielectric cavities, which can be called 'regular' cases 
since their closed counterparts (billiard problems) are in- 
tegrable. To our knowledge no analytical solution has 
been proposed so far for these cavity shapes in the open 
case. Nevertheless, their dielectric spectrum shows some 
characteristic features specific to integrable systems. 

The only example of 2D integrable dielectric cavity is 
the circular one (see e.g. [22]). Therefore its resonances 
are organized in regular branches labeled by well defined 
quantum numbers. For the above mentioned 'regular' 
cavities with relatively small refractive index it appears 
that the resonances still follow similar branch structure 
(see Figs. [3j [6j [9j . This is surprising as these dielectric 
problems are not integrable and strictly speaking there 
is no conserved quantum number. This unusual regu- 
larity can be described by the superscar approximation 
proposed in [To| • The detailed discussion of the modified 
superscar model and its application to these problems 
will be given elsewhere (23| . 




FIG. 2: (a) Two representations of the same diamond periodic 
orbit. Definition of the incident angle %. (b) The dashed part 
corresponds to the fundamental domain used for numerical 
simulations. The Fabry-Perot and diamond periodic orbits 
are drawn, as well as their restriction to the fundamental do- 
main (dotted lines). 



A. The square cavity 

The simplest example of regular cavities is the dielec- 
tric square where the inside billiard problem is straight- 
forward. For instance, for Dirichlet boundary conditions, 
the eigenenergies and eigenfunctions are the following: 

7T 2 

k lm = ~2(P 2 +rn 2 ) , 

2 . fp7rx\ . /miry\ 
ipp,m{x,y) = -sin sin , 

where a is the side length and p and m are two positive 
integers. The outer scattering problem even with the 
Dirichlet boundary conditions is more difficult as it cor- 
responds to a pseudo-integrable problem (see e.g. [24|) 
and no explicit analytical solution exists [25j. 



1. Numerics 

The solutions of the square dielectric problem can be 
divided into four symmetry classes corresponding to wave 
functions odd (— ) or even (+) with respect to the diago- 
nals y = x and y = — x. For instance, the notation ( — h) 
means that the wave function is odd with respect to the 
diagonal y = x and even with respect to the other. Each 
symmetry class reduces to a quarter of a square (dashed 
part in Fig. (2}d) with Dirichlet, (— ), or Neumann, (+), 
boundary conditions along the diagonals. The ( — +) and 
(H — ) symmetry classes are equivalent. Fig. [3] shows the 
resonance spectrum for all symmetry classes and n = 1.5. 
Fig. H] displays some typical quasi-stationary states for 
the ( ) symmetry class from different parts of the spec- 
trum. Notice highly unusual regularity of the spectrum 
and wave functions for this shape. 

As for the case of the circular cavity (see e.g. (22|), 
the imaginary part of the dielectric square resonances is 
bounded by the losses j m ax of the periodic orbit with the 
highest losses, which is here the Fabry-Perot indicated in 
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FIG. 3: (Color online) Resonance spectrum for the dielec- 
tric square with n — 1.5. The crosses, squares, and circles 

correspond respectively to ( — h), (++), and ( ) symmetry 

classes. The position of the horizontal dashed line is given by 

G3. 



Fig. [3d with 
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-1.073 



=1.5 



(13) 



The counting function N(k) gives the mean number of 
resonances with a real part less than k in the strip defined 
by (jU). The Weyl-type formula (|6j) estimates its growth 
when k — >• oo. We checked this prediction for different 
values of the refractive index and each symmetry class. 
For n = 1.5 the results of the numerical fit to the data 
computed from the spectrum in Fig. [3] are the following 



(— ) Nfu(k) 
(++) N fit (k) 
(-+) N fit (k) 



n 2 



16tt 



(kay 



0.0866 ka- 0.0456, 
0.2489 ka- 1.737, 
0.0806 ka- 2.012. 



Here we fixed the coefficient of the quadratic term and 
fitted the linear and constant terms from the numerical 
data. 

The predictions of ([7]) which take into account the 
Dirichlet or Neumann boundary conditions [36| on two 
of the sides of the fundamental domain are given by the 
following expressions calculated for n = 1.5 
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(-+) a th = 
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0.0816. 



n=1.5 




(c) 



FIG. 4: Wave functions for the dielectric square with n = 1.5 

and ( ) symmetry class, (a), (b), and (c) correspond to 

ka = 98.36 - 0.0071 i, ka = 98.78 - 1.061 i, and ka = 98.96 - 
0.998 i. Greyscale, black representing maximal values of \ip\ 2 . 



The predicted values are in good agreement with our nu- 
merical calculations. More precisely, each residue (the 
difference between N(k) and Nf it (k)) oscillates around 
zero as evidenced in Fig. [5] To ensure the efficiency of 
formula ([6]) we calculated the spectrum also for n = 1.35 

with ( ) symmetry class, see Fig [6] The quadratic fit 

for the counting function gives now 



Nfu(k) 



16tt 

to be compared with: 

r(n) — n\/2 



(ka) 2 - 0.0709 ka - 2.052 



a 



th 



47T 



-0.0713. 



= 1.35 



The residue stays also close to zero (see insert in Fig. [6]). 
In all investigated examples, the agreement between pre- 
diction and numerics is better than 2% for the linear 
term. 
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FIG. 5: Residue N(k) — Nf it (k) for the dielectric square (a — 

1) with n — 1.5 and different symmetry classes : (a) ( ), 

(b) (++), and (c) (-+). 
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FIG. 6: Resonance spectrum for the dielectric square with 

n — 1.35 and ( ) symmetry class. The positions of the 

horizontal dashed lines are given by ([T3]) . Insert: N(k) — 
Nfit(k) with a = 1. 



The oscillatory part of the trace formula is checked 
as well. In the square, periodic orbits form families. 
Thus their weighting is predicted by (fTT]h which im- 
plies that the spectrum is dominated by the diamond 
periodic orbit (see Fig. [2^). Actually the weighting co- 
efficient of this orbit is calculated as follows: it covers 
the whole cavity (Aiiamond = & 2 ), its length is short 
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FIG. 7: Comparison of \d(l)\ 2 calculated from the numerical 

data of the dielectric square for ( ) symmetry class and 

n — 1.5 (continuous line), and predictions from trace formula 
(crosses) considering the m th repetitions of half the diamond 
periodic orbit (/ = may/2) with a coefficient |<i(/)| 2 oc 1/ra. 



(^diamond = 2ay/2), and for n > \/2 there is no refrac- 
tive loss ( | Rdiamond | = 1) in the deep semiclassical limit 
k — >• oo. For illustration, it is worth comparing with 
the c p coefficient of the Fabry-Perot periodic orbit (37| : 
^4 FP = a 2 , £ FP = 2a, and it^p = [(n - l)/(n + l)] 2 . 
Then for n = 1.5, \c FP /c diarnond \ ~ 0.05 mainly due to 
the prominent influence of refractive losses. 

The agreement with numerical simulations is checked 
via the length density of the dielectric square which is 
computed from the numerical spectrum using ([T2]) . The 
results are plotted in Fig. [71 The length density is highly 
peaked at I = a\[2 and at its harmonics (l m = may/2, 
with m G N*). Only half of the diamond orbit length 
appears, since the length density is calculated for a sin- 
gle symmetry class and Fig. [2]o shows that the diamond 
periodic orbit is twice shorter if restricted to the dashed 
area. If the length density had been performed with the 
four symmetry classes, it would have been peaked at the 
full diamond length (2ay/2) and at its harmonics. The 
same appears with experiments as shown below. 

The agreement between numerics and predictions from 
trace formula JTT} is quite good as well when comparing 
the ratios of the harmonics. Actually these harmonics 
can be identified as repetitions of the diamond periodic 
orbit (C = m£diamond) and thus formula ([TT]) predicts 
that the |c p | 2 should decrease like 1/ra. This prediction 
is shown by crosses in Fig. [71 From numerics, we received 
harmonics a little bit smaller than predicted which is 
natural as the Fresnel reflection coefficient ([TO]) does not 
take into account correctly a leakage through a dielectric 
interface of finite length. 



2. Experiments 

The prevalence of the diamond periodic orbit was al- 
ready experimentally demonstrated with organic micro- 
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FIG. 8: (Color online) (a) Experimental spectra from a square 
micro-laser with a — 120 /mi. Top: a — 45°, bottom: a — 
0°. (b) Fourier transform of the spectra in (a). The arrow 
indicates the predicted optical length of the diamond periodic 
orbit: 2\^2an. 
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FIG. 9: Resonance spectrum of the dielectric rectangle with 
p — 2, (++) symmetry class, and n=1.5. The positions of the 
horizontal lines are given by (|14|) and (|15j) . 



diamond length) is slightly higher than the first one. This 
is due to the presence of a residual comb visible in the 
spectrum (Fig. [8^ top). 

In this section, we have shown that the spectral prop- 
erties of the dielectric square (density of states, resonance 
losses, laser spectra) are controlled in a first approxima- 
tion by classical features and symmetry classes. For low 
refractive indices (what we studied), the diamond peri- 
odic orbit plays a prominent role. 



B. The dielectric rectangle 



lasers [loj and micro- wave cavities [26j. Here we would 
like to stress that sometimes the experimental spectrum 
reveals half the diamond periodic orbit instead of the 
full one due to a selection of symmetry classes. This 
phenomenon is illustrated in Fig. using the pump po- 
larization as a control parameter. Actually the DCM 
molecule (the laser dye) is more or less 'rod-like' and thus 
conserves (in a way) the memory of the pump polariza- 
tion which can be monitored at will without modifying 
other parameters. A study of the pump polarization in- 
fluence will be published elsewhere [l9| . 

The Fourier transforms of the spectra in Fig. [8^ are 
plotted in Fig. [8]o. Let us call a the angle between the 
pump polarization (which lies in the plane of the layer) 
and the direction of observation. For a = 0°, the first 
harmonic of the Fourier transform is peaked at the di- 
amond optical length. But for a = 45°, only one peak 
out of two appears in the spectrum, therefore the Fourier 
transform peaks at half the diamond optical length. It 
could be noted that the second harmonic (at the actual 



We repeat the same steps as in the previous Section 
but for a rectangular cavity so as to monitor eventual 
changes when breaking the square symmetry. Let us call 
p = L/l the ratio between the larger and smaller sides. 
We will focus on the case p = 2. 

1 . Numerics 

We restrict ourselves to the (++) symmetry class with 
respect to the perpendicular bisectors of the sides (see 
dashed area in Fig. HTb). Fig. [9] shows the resonance 
spectrum for this case. The lower bound 7 max of the 
imaginary parts of the resonances is related to the life- 
time of the classical orbit bouncing off perpendicularly 
the longest side of the rectangle: 

Imax = "In ( ^— ^ ) (14) 

Another horizontal line is plotted in Fig. [9j which corre- 
sponds to the lifetime of the orbit bouncing off perpen- 
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FIG. 10: Residue N(k) — Nfu(k) for the dielectric rectangle 
with p = 2, I = 1, (++) symmetry class, and n— 1.5. 





■ 








I = geometrical length 

(a) 

FIG. 11: (a) Length density for a dielectric rectangle calcu- 
lated from the spectrum in Fig. [9] The expected position 
of the double diamond (resp. the stretched diamond) is in- 
dicated by a cross (resp. a triangle), (b) Main periodic or- 
bits contributing to the resonance spectrum: double diamond 
(dashed line) and stretched diamond (continuous line). The 
dashed area corresponds to the fundamental domain. 



dicularly the smallest side of the rectangle: 



72 = — In 

pn 



(15) 



As for the square cavity, prediction (j6j) is checked. The 
best quadratic fit of the counting function N(k) com- 
puted numerically from data of Fig. [9] is (with n = 1.5): 

N fit (k) = n 2 — (klf + 0.2894 kl - 3.7346 . (16) 

167T 

The prediction for the linear term is: 
r(n) 



th 



8tt 



n=1.5,/9=2 



: 0.3014 



(17) 



which shows a good agreement. The difference between 
the numerically computed N(k) and its best quadratic 
fit ([II]) is shown in Fig. [HB 

The length density d(l) defined by ([T2]) is shown in 
Fig. [TTa , and is peaked at the lengths of the double di- 
amond orbit and the stretched diamond orbit (both dis- 
played in FigfTTb). 

Eventually some wave functions associated to reso- 
nances from different parts of the spectrum are shown 
in Fig [H 
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FIG. 12: Quasi- stationary states for the dielectric rectangle, 
n = 1.5, (++) symmetry class, (a) kl = 70.04 - .038 i, (b) 
kl = 70.26 - .33 i. (c) kl = 70.45 - .92 i. Greyscale, black 
representing maximal values of \ip\ 2 . 



2. Experiments 



Fig. [T3k presents a typical experimental spectrum from 
a rectangular microlaser with p = 2. Its Fourier trans- 
form plotted in Fig. [T3b is peaked at the length of the 
double diamond periodic orbit (see insert), in agreement 
with numerics and predictions. This experimental ob- 
servation is very robust whatever the parameter being 
used: direction of emission, pump intensity, and pump 
polarization. For illustration, a comparison between the 
measured and expected optical lengths is presented in 
insert of Fig. [T3b for various cavity sizes. For complete- 
ness, it should be noticed that the Fabry-Perot along the 
longest axis appears if observed in its specific direction 
and pumped with a favorable polarization. 
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FIG. 13: (Color online) (a) Experimental spectrum of a rect- 
angle microlaser with p — 2 and / = 70 pm. (b) Fourier trans- 
form of the spectrum in (a). Insert: Comparison between the 
optical lengths inferred from experiments and those expected 
for the double-diamond for different cavity sizes. 
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FIG. 14: Resonance spectrum for the dielectric ellipse with 

p — 2, ( ) symmetry class, and n = 1.5. The position of 

the horizontal dashed line is given by JT5J). 



C. The dielectric ellipse 




The ellipse can also be considered as a 'regular shape', 
since the interior billiard problem is separable [27]. Let 
us call b (resp. a) half the length of the minor (resp. ma- 
jor) axis. Here we consider only the ratio p = a/b = 2, 
however the computations for other p values give simi- 
lar results. Here we will restrict ourselves to the ( ) 

symmetry class, i.e. the function vanishes along both 
symmetry axis of the ellipse. 



1 . Numerics 

Fig. HH shows the resonance spectrum for p = 2. As for 
the square cavity, it looks quite regular while the problem 
is not separable. Similarly, the imaginary parts of the 
resonances are bounded by the losses of the Fabry-Perot 
periodic orbit (along the minor axis): 



2n 



In 



n - 



1 



(18) 



FIG. 15: N(k)-N fit (k) for the dielectric ellipse (b = 1) with 
p = 2 and ( ) symmetry class. 



Formula (|6]) for the counting function is checked with 
the same protocol as before. The residue between nu- 
merics and fit is plotted in Fig. [15] and oscillates around 
zero. Moreover the linear term of the regression 

a fit = -0.0866 

agrees well with the prediction: 

th _ r(n)pE(e) - n(l + p) 



a 



47T 



n=1.5,/9=2 



-0.0872 



where E(z) is the complete elliptic integral: 

r ir/2 



E(z) = / yj\-z 2 sin(£) 2 d£ , 
Jo 



(19) 



and e = yl — (b/a) 2 is the eccentricity of the ellipse. 
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FIG. 16: Wave functions for the ellipse with p = 2, n = 1.5 

and ( ) symmetry, (a) Whispering gallery mode, kb = 

22.00 - 0.02 i. (b) Bouncing ball mode, kb = 22.30 - 0.36 i. 
Greyscale, black representing maximal values of l^l 2 - 
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FIG. 17: Fourier transform of the spectrum in insert. The 
arrows indicate the predicted optical length of periodic or- 
bits: single arrow for the Fabry-Perot along th e major axis 
4a n and double arrow for the rectangle Abn^/l + p 2 . In- 
sert: Experimental spectrum from an elliptical micro- laser 
with p = a/b = 2 and b = 50 pm registered in the direction 
parallel to the major axis. 



The oscillatory part will be postponed to a future pub- 
lication. We already note that the wave functions display 
in general two kinds of behavior for the wave inside the 
cavity: either "whispering gallery modes" or "bouncing 
ball modes" - like as is (rigorously) the case for the el- 
liptic billiard. Fig. [16] presents examples of such wave 
functions, which correspond to resonances from different 
parts of the spectrum. 



2. Experiments 

Experiments provide similar insights onto the domi- 
nant resonance features. The insert in Fig. [17] shows 
a typical spectrum from an elliptical micro-laser with 
p = 2, while its Fourier transform is plotted in the main 
window. Its first harmonics presents two main peaks with 
positions corresponding quite well to the optical lengths 
of two periodic orbits: the rectangle and the Fabry-Perot 
along the major axis. The deviation is less than 3 %, 
which is the experimental inaccuracy. It should be noted 



that for p = 2 the length of the Fabry-Perot along the 
major axis is equivalent to the second repetition of the 
Fabry-Perot along the minor axis, also called bouncing 
ball. 



IV. PSEUDO-INTEGRABLE SYSTEM: THE 
DIELECTRIC PENTAGON 

Similar studies were performed for the dielectric pen- 
tagon and hexagon, which are particularly interesting 
systems since they contain diffracting angles: mir/n with 
m, n co-prime and m > 1. Following Richens and Berry 
|24{ these systems are called pseudo-integrable, since 
their classical flow is confined to a surface as for inte- 
grable systems but because its genus is bigger than 1 
they cannot be classically integrable. 

Here we only consider the dielectric pentagon, though 
every conclusion also applies to the hexagon mutatis mu- 
tandis. Below we present the results for the ( ) sym- 
metry class, which means that the associated wave func- 
tions vanish along each symmetry axis of the polygon. R 
stands for the radius of the outer circle of the pentagon 
and a = 2i?sin(7r/5) for its side length. 



A. Numerics 

Fig. [T8a shows the resonance spectrum of a dielectric 

pentagon for ( ) symmetry and n = 1.5. Again the 

imaginary part of the resonances is bounded and this 
lower bound jrnax can be estimated from the refractive 
losses of the periodic orbit drawn in Fig. [T8b and d which 
presents the highest losses (shortest lifetime): 



2, n / 7T \ 
Imax = —7 m RtM I I 

TIL 



-0.630 



(20) 



VlO/ 

where I = 2R (sin(27r/5) + sin(7r/5)). 

Several wave functions are displayed in Fig. [19] It is 
important to stress the existence of different types of res- 
onances. The ones with low losses are related to 'whis- 
pering gallery-like' modes, see Fig. [T9a. In [10], this ob- 
servation was used to build a superscar approximation 
of these resonances. Notice the more complex pattern 
in (b) for this wave function corresponding to kR with 
rather large imaginary part, and the scarring by the orbit 
in Fig. [TSll for the wave function in Fig. [T9b . 

The Weyl law ([6]) is checked as above by fitting the 
counting function N(k) which gives: 



Nfu(k) 



n 2 sin(27r/5) 

167T 



(kR) 2 -0.1681 hR+ 0.1077 , 



in good agreement with the prediction for the linear term: 
th r(n) sin(7r/5) — n [1 + cos(7r/5)] 



47T 



-0.1680 



n=1.5 



Moreover the residue N(k) — Nf it (k) oscillates around 
zero as expected (see Fig. [T8b). 
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FIG. 18: (a) Resonance spectrum for the dielectric pentagon 
with ( ) symmetry and n — 1.5. The position of the hor- 
izontal dashed line is given by (|20p . (b) N(k) — Nfu(k) cal- 
culated from the numerical spectrum in (a) with R = 1. (c) 
A part of the periodic orbit with the shortest lifetime, (d) A 
representation of the periodic orbit in (c). 



The length density plotted in Fig. [20b evidences that a 
few periodic orbits mostly contributes to the oscillatory 
part of the trace formula. First the density of orbit length 
is peaked at a length corresponding to the 'double pen- 
tagon' periodic orbit, which is depicted in Fig.l20k and b. 
This orbit is confined by total internal reflection for our 
value of the refractive index and lives in family contrary 
to the isolated 'single pentagon' orbit. The length of the 
single pentagon periodic orbit, once folded in the fun- 
damental domain shown in Fig. [20k . is / = i?sin(27r/5). 
The length of the double pentagon periodic orbit is twice 
longer. The vertical lines in Fig.l20b indicate the theoret- 
ical lengths of the m th repetition of the single pentagon 
periodic orbit: l m =ml. As expected, the length density 
is mostly peaked at the l m positions with even m. 

Second, it is worth noting that the amplitude of the 
peaks does not clearly decay as for the square in Fig. [71 
For m = 7 for instance, the peak amplitude is unexpect- 





FIG. 19: Wave functions of the dielectric pentagon with ( ) 

symmetry and n = 1.5. (a) kR = 91.35 - 0.038 i, (b) kR = 
89.06 - 0.453 i, and (c) kR = 92.56 - 0.628 i. Greyscale, black 
representing maximal values of l^l 2 - 



edly high for a repetition of a given orbit. This happens 
when the length l m is close to the length of a diffracting 
orbit. Here, the repetitions m = 7 and m = 9 are in fact 
of similar lengths than the orbits illustrated in Fig. [2TJ 
The treatment of such diffracting correction (see a simi- 
lar discussion in [28| for billiards) is beyond the scope of 
this paper as it requires the local exact solution of the 
diffracted field by a dielectric wedge. 



B. Experiments 

A typical experimental spectrum from a pentagonal 
micro-laser is plotted in Fig. [22k . Its Fourier transform 
presented in Fig. [22b is mostly peaked at the length of the 
double pentagon periodic orbit. The single pentagon is 
visible as well, which can be directly noticed on the spec- 
trum made of two combs of different amplitudes. In |10l |. 
we reported an experimental spectrum from a pentago- 
nal micro-laser where both combs had similar amplitudes 
and therefore its Fourier transform did not present any 
peak at the length of the single pentagon. The parame- 
ters which control the relative amplitudes of the combs 
have not been identified yet, but it is clear that the etch- 
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FIG. 20: (Color online) (a) 'single pentagon' periodic orbit 
in solid line and 'double pentagon' periodic orbit in dotted 
line. The dashed area indicates the fundamental domain used 
for simulations, (b) A representation of the double pentagon 
periodic orbit, (c) Length density for the dielectric pentagon 
calculated from the numerical data plotted in Fig. [T8h . The 
vertical lines in the bottom part of the graph indicate the 
repetitions of the single pentagon orbit. 
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Diffracting orbits in the pentagonal cavity. 
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FIG. 21 

(a) l a = R a /130 + 22 y^ /2 

(b) l b = RV210 + 38^5/2, 9l/l b 
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FIG. 22: (a) Experimental spectrum from a pentagonal micro- 
laser with a = 90 /im. (b) Fourier transform of the spectrum 
in (a). The single (resp. double) arrow indicates the expected 
position of the single (resp. double) pentagon periodic orbit. 



A. Numerics 

For simplicity we only consider the ( ) symmetry 

class, which means that the associated wave functions 
vanish along both symmetry axis of the stadium. The 
resonance spectrum for p = 1 is shown in Fig. [23k . As for 
the other cavities, the imaginary part of the resonances is 
bounded and its lower bound jmax can be estimated from 
the refractive losses of the periodic orbit which presents 
the highest losses, i.e. the Fabry-Perot along the small 
axis (also called bouncing ball orbit): 



-_Li ( n ~ l 

lmax ~ 2n U+l 



(21) 



ing quality is a key point. 



Again, some wave-functions are presented in Fig. EU 
representative of different parts of the spectrum. 

The counting function N(k) was computed from the 
numerical data shown in Fig. [23k and the best fit gives: 



V. CHAOTIC DIELECTRIC CAVITIES 

Eventually we applied the same ideas to an archetypal 
chaotic cavity, the Bunimovich stadium, which is made 
of a rectangle between two half circles (see Fig. [23b for 
notations), and investigated various deformation defined 
by the parameter p = L/R. 



Nfu(k) = — [p + - J (kR) 2 - 0.145 kR - 4.042 , (22) 

to be compared with the prediction ([7]) for the linear 
term: 



th 



1 

47T 



-n(2 + p) 



n—1.5,p=l 



.148 . 



(23) 
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FIG. 23: (a) Resonance spectrum for the dielectric stadium 

with p — 1, ( ) symmetry class, and n — 1.5. The position 

of the horizontal dashed line is given by ([2T]) . (b) Notations 
for the Bunimovich stadium, (c) N(k) — Nf it (k) calculated 
from the numerical simulations in (a) with R = 1. 



The oscillatory part of the trace formula is also 
checked, plotting the length densities calculated from nu- 
merical spectra for several shape ratios p. The curves pre- 
sented in Fig. [25] are peaked at different positions which 
could be assigned to periodic orbits. To predict which pe- 
riodic orbits should mainly contribute to the length den- 
sity, we calculated their weighting coefficient from for- 
mula (j9j). The considered orbits are drawn in Fig. [26j 
their geometrical length is plotted in Fig. [27J and their 
coefficient (amplitude) versus p in Fig. [28k . Note that 
orbits 5 and 6 obey geometrical constraints such as they 
do not exist for p < 1. 

The length densities are calculated from the numerical 
spectra on a finite number of resonances (finite range of 
Re(kR)), thus some finite size effects do play a role and 
must be taken into account evaluating the amplitudes of 
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FIG. 24: Wave functions of the dielectric stadium with p = 1, 

( ) symmetry class, and n — 1.5. (a) kR = 34.48 — 0.11 i. 

(b) kR = 34.55 - 0.31 i. (c) kR = 34.48 - 0.45 i. Greyscale, 
black representing maximal values of \i/j\ 2 . 



the periodic orbits. One of the main effect here comes 
from the curvature correction in the reflection coefficient. 
If the dielectric boundary is curved enough compared to 
the wavelength then there is a quite important correc- 
tion to the standard Fresnel coefficients [29] . We use the 
following formula to take into account the curvature cor- 
rection: 



o(c) _ 

n TM — 



H 



m 



(i)' 



(24) 



where x = kR and m = nx cos x- The Fresnel coefficient 
for a straight interface (fTU|) is recovered for large kR, see 
Fig.[29l In Fig. [28b, the amplitudes of the periodic orbits 
are plotted taking into account this curvature correction. 
Notice the important differences with Fig. [28k . 

Using Fig. [28b we can give a quantitative estimate of 
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/ = geometrical length 

(c) 

FIG. 25: Length densities calculated from numerical spectra 

for dielectric stadiums with ( ) symmetry class and n — 

1.5. (a) p = 0.3. (b) p = 0.75. (c) p = 1.25. The vertical 
dotted lines indicate the expected positions of periodic orbits 
(numbering according to Fig. I26[) . 




FIG. 26: Some periodic orbits of the stadium. 




p 

FIG. 27: Geometrical length of the periodic orbits plotted in 
Fig. 1261 Numbering according to Fig. [26] 
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FIG. 28: Amplitude of the orbits listed in Tab. [26] for the 

( ) symmetry class following formula J9| (a) Semiclassical 

regime (i.e. no size effect taken into account), (b) Idem with 
curvature correction calculated for kR — 25. 
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X 

FIG. 29: Modulus of the TM reflection coefficient as a func- 
tion of the incidence angle (in degrees). Full line: Fresnel 
coefficient (|1Q|) for a straight boundary. Dotted lines: (|24|) for 
kR 100, 50, 25. 



the periodic orbits which mostly contribute to the length 
density. For p = 0.3 (Fig. [25k). the length density is 
peaked around the orbits 3 (Z ~ 3.28) and 1 (Z ~ 3.43), 
and orbits 2 and 4 with respective length Z ~ 3.6 and 
Z ~ 3.54 cannot be separated. For p = 0.75 (Fig. [25b), 
the orbits 1 and 4 with respective length Z ~ 4.33 and 
Z c± 4.40 interfere. The line at Z 4.5 stands for orbit 2. 
Eventually for p — 1.25 (Fig. [25b), the two orbits 1 and 
4 (Z ~ 5.33 and Z ~ 5.38) interfere. Again we drew a line 
for orbit 2 at Z ~ 5.5. A peak can also be seen for orbit 
6 (Z ~ 6.04). From these examples it appears that the 
agreement between theory and numerics is qualitatively 
good. 



B. Experiments 

Spectra were recorded for the shape ratios used in 
numerical simulations and for various R (see [30| for 
p = 0.5). A typical experimental spectrum from an or- 
ganic micro-stadium is plotted in Fig. [30] and typical 
Fourier transforms in Fig. [3TJ They look similar to those 
of the shapes studied in the previous sections. However 
the main difficulty to face when studying stadiums is the 
large number of orbits with close lengths. Would it be 
numerically or experimentally, it is thus difficult to as- 
sign peaks in the length density. Then, to check formula 
© against experiments, it was decided to compare the 
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FIG. 30: Experimental spectrum of a stadium micro-laser 
with p — 0.75 and R = 50 pm. Its Fourier transform is plotted 
in Fig. GEO). 
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FIG. 31: Normalized Fourier transform of experimental spec- 
tra from stadium micro-lasers with (a) R = 80 /im, (b) 
R = 50 /im, and (c) R = 45 pm. The vertical line corre- 
sponds to the position of the perimeter. The arrows indicate 
the expected positions of some orbits, for comparison: (a) or- 
bits 1 and 2, (b) orbit 1, and (c) orbits 1 and 5 (numbering 
according to Fig. 126)) . 

positions of the peaks to the length of the perimeter. Ac- 
tually some crossing orbits like orbits 5 and 6 are longer 
than the perimeter and, due to geometrical constraints, 
do not exist for p < 1. Moreover according to Fig. [28k , 



which corresponds to the semi-classical limit and so to ex- 
perimental conditions, their amplitudes (0 are the high- 
est when they appear. So we expect Fourier transforms 
peaked at positions shorter than the perimeter for p < 1 
and longer for p > 1, and this is evidenced in Fig. [3TJ 



VI. CONCLUSION 

In this paper we have shown numerical and experi- 
mental results concerning the trace formula for dielectric 
cavities. For convex cavities and TM polarization the res- 
onance spectrum can be divided into two subsets. One 
of them, the Feschbach (inner) resonances which are rel- 
evant for experiments, is statistically well described by 
classical features: the periodic orbit with the shortest 
lifetime for the lower bound of the wave number imag- 
inary parts, the Weyl's law for the counting function, 
and the weighting coefficients of periodic orbits for the 
length density. The formulae we derived, based on stan- 
dard expressions used in quantum chaos and adapted to 
dielectric resonators, give an accurate description of these 
spectral properties. 

Those formulae have been checked for various resonator 
shapes. For "regular shape" (i.e. the analogous billiard 
problem is separable) and for small index of refraction, 
the resonance spectrum presents a branch structure as 
if the dielectric problem was separable. When the cor- 
responding billiard problem is not integrable, the usual 
Weyl's law still occurs. The oscillating part can be also 
explained by taking into account the shortest periodic 
orbits. In the chaotic case the correspondence is how- 
ever more difficult to claim quantitatively as finite size 
effects play a quite important role because of periodic 
orbits with close lengths. 

The main result of the paper is the demonstration that 
dielectric cavities widely used in optics and photonics can 
be well described using generalizations of techniques from 
quantum chaos. 

This study raises many open problems and we would 
like to mention some of them. The formulae checked here 
gave accurate predictions for TM polarization. For TE 
polarization, due to the existence of the Brewster angle 
where the reflection coefficient vanishes [3l| , the situation 
is less clear and requires further investigations. The next 
step should be to treat carefully the diffraction on dielec- 
tric wedges, which is still an open problem [32|. A related 
question is to improve the accuracy of the standard effec- 
tive index theory, since the separation into TE and TM 
polarizations is precisely based on this 2D approximation. 
Experimental data (i.e. real 3D systems) indeed revealed 
departures from this model [H, |26[ . These questions are 
related to the wave functions (see e.g. [loj ) and far field 
patterns, which are also of great interest, especially for 
applications. 
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FIG. 32: Plot of r(n). 
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Appendix A: Weyl's law 

The Appendix deals with the derivation of formula (j6j) 
using an alternative method than in [3[. Start from the 
definition of the Green function G: 



(A* + n(x fk 2 ) G(x , y) = 5{x - y) 



(Al) 



where n(x ) is equal to n (resp. 1) when x is inside (resp. 
outside) the dielectric cavity. Moreover, for TM modes, 
G(x,y) and its normal derivative are continuous along 
the boundary of the domain. Taking the trace of it gives 
the density of states through the Krein formula (see Q 
and references therein): 

daii(E) - do(E) = (A2) 
— — j Im [n(x ) 2 G(x, x ) — Gq(x : x )] dx 

where the integral runs over the whole 2D plane. do(E) is 
the density of states of the free space and Go(x,y) stands 
for the Green function of a free particle in the plane: 



G (x,y) = ^H^\k\x-y\) 
1 / 4+oc e ip ( Xl ~ yi ^> 



(A3) 



47ri 



Vk 2 



>i\A 2 -P 2 \ X2 ~y 2 \dp (A4) 



with x — {x\,x<2) and y — (2/1,2/2)- H^\z) is the Han- 
kel function of the first kind and Eq. (|A4|) implicitly as- 
sumes that k has a small positive imaginary part. It is 
worth noting that d a u{E) is the Krein spectral shift func- 
tion which is different from the spectral density d(E) dis- 
cussed in Eq. (J5j). It will be explained how to get it at 
the end of this Section. 



1. Derivation of the first two terms 

The leading term of the Weyl's law is obtained when 
substituting G to G in Eq. (IM1) and using (|X3j) : 

d all (E)-d (E)^(n 2 -l)^- (A5) 

47T 

where A is the area of the domain filled with the dielectric 
material. 

The first remaining term in (j6j) comes from the pres- 
ence of the boundary. Thus it is first necessary to solve 
the elementary problem of a plane wave reflecting on an 
infinite straight dielectric boundary, which can be derived 
through standard methods. Then, using expression (jA4p , 
the Green function for both x and y inside the dielectric 
can be written: 



G(x,y) = 



1 



+ 00 e ip(x 1 -y 1 ) r —- 

.\ P W n k p F2-S/2I 



47riJ_ 00 y/n 2 k 2 - p 2 1 " 

R{p) e W" 2 k 2 -p 2 |x2+?/2|] dp 



(A6) 



where p is the tangential component of the momentum 
and 



R(p) = 



\J n 2 k 2 — p 2 — \/k 2 — p 2 



^Jn 2 k 2 — p 2 + ^/k 2 — p 2 



(A7) 



Similarly one gets the Green function when the argu- 
ments are outside the dielectric: 



G{x,y) 



1 

47ri 



+00 ip(x 1 -y 1 ) 



> \]k 2 - p 2 



D i^k 2 -p 2 \x 2 -V2\ 



(A8) 



As usual the trace of G is computed using local coordi- 
nates. The surface term (|A5|) is recovered from the first 
terms of ()A6|) and (|A8|h so we focus now on their second 
terms only. The integration along the boundary gives 
the length factor C. For the transverse coordinate, the 
boundary is approximated locally by its tangent plane, 
and then (|A6|) and (|A8|) are used. After the convenient 
Wick rotation p — >• — itk, the boundary contribution from 
inside is: 



8tt 2 /c 



/+00 
-00 



R(-itk) 
n 2 +t 2 



dt . 



Similarly from (|A8|I the boundary contribution is: 

r+0 ° R(-itk) 



(out) 



1 r°° r(_ 



-t 2 



-dt . 



(A9) 



(A10) 



Putting together (|A9|) and (|A10p back to (|A2|) . one gets 
the first two terms of the Weyl expansion: 



dall{ E)-d (E)c(n 2 -l)£+ {f{n) 1)£ 



47T 



Snk 



, (All) 
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where, noting R(t) instead of R(—itk): 

,2 ,+oo df 1 ,+00 df 



f(n) = 1 + 



t 2 + n 2 y ' ir 



t 2 + l 



r(n) is plotted in Fig. 



R(t) . 
(A12) 



From the density of states to the count of 
Feschbach resonances 



As mentioned above, the quantity entering the Krein 
formula is not exactly the spectral density (j5j), but they 
can be related heuristically. The spectral shift function, 
d a u(E) in (jA2|h is related with the determinant of the full 
S'-matrix for the scattering on a cavity, while d(E) in (j5j) 
is the spectral density of Feschbach (inner) resonances 
which are poles of this ^-matrix. In addition to these 
poles, the determinant of the S'-matrix may have poles 
associated with shape resonances (which we do not take 
into account) and then an additional phase factor d s (E): 



It is natural to assume (and can be checked for dielec- 
tric disk) that for all outside structures the corresponding 
wave functions are almost zero inside the cavity and on 
its boundary Therefore the S- matrix phase, d s (E), as- 
sociated with such functions in the leading order is the 
same as for the outside scattering on the same cavity 
but with the Dirichlet boundary conditions. The Weyl 
expansion for such scattering is known (see 0, Q and 
references therein): 



ds{ E)- M E)^-±-^ k . (A14) 



Subtracting (|A14|) to (|A11|) gives the desired result for 
the spectral density, taking into account only Feschbach 
resonances: 



d(k) = 2kd(E)~n 2M ' f(n) 



2tt 



47T 



C . 



(A15) 



d all {E) = d{E) + d s {E). 



(A13) Then Eq. (j6]) is recovered by integration. 
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